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- - - ^ We discuss the dynamical system approach applied to Higher Order Theories of Gravity. We 

, show that once the theory of gravity has been specified, the cosmological equations can be written 

3C ' as a first-order autonomous system and we give several examples which illustrate the utility of our 

I method. We also discuss a number of results which have appeared recently in the literature. 

^ ; I. INTRODUCTION 

Although there are many good reasons to consider General Relativity (GR) as the best theory for the gravitational 
interaction, in the last few decades the advent of precision cosmology tests appears more and more to suggest that 
this theory may be incomplete. In fact, besides the well known problems of GR in explaining the astrophysical 
phenomenology (i.e., the galactic rotation curves and small scale structure formation), cosmological data indicates 
an underlying cosmic acceleration of the Universe which cannot be recast in the framework of GR without resorting 
^ : to additional exotic matter components. Several models have been proposed PJ in order to address this problem 
and currently the one which best fits all available observations (Supernovae la Cosmic Microwave Background 
anisotropies [1], Large Scale Structure formation [J, baryon oscillations 0], weak lensing [6j]), turns out to be the 
Concordance Model in which a tiny cosmological constant is present 0] and ordinary matter is dominated by a Cold 
Dark component. However, given that the A-CDM model is affected by significant fine-tuning problems related to 
the vacuum energy scale, it seems desirable to investigate other viable theoretical schemes. 

It is for these reasons that in recent years many attempts have been made to generalize standard Einstein gravity. 
Among these models the so-called Extended Theory of Gravitation (ETC) and, in particular, non-linear gravity 
theories or higher-order theories of gravity (HTG) have provided interesting results on both cosmological [sl. [ol liol [Til. 
H, cind astrophysical [ll|, [U scales. These models are based on gravitational actions which are non-linear in the 
' Ricci curvature R and/ or contain terms involving combinations of derivatives of R [H, [l^, [l3|- The peculiarity of 
I these models is related to the fact that the gravitational field equations can be recast in such a way that the higher 
order corrections provide an energy - momentum tensor of geometrical origin describing an "effective" source term on 
the right hand side of the standard Einstein field equations 8, ^10!|. In this scenario, the cosmic acceleration can be 
shown to result from such a new geometrical contribution to the cosmic energy density budget, due to higher order 
^ corrections to the Hilbert-Einstein Lagrangian. 

Because the field equations resulting from HTG are extremely complicated, the theory of dynamical systems provides 
a powerful scheme for investigating the physical behaviour of such theories (see for example |l8l.[l9j). In fact, studying 
cosmologies using the dynamical systems approach has the advantage of providing a relatively simple method for 
obtaining exact solutions (even if these only represent the asymptotic behavior) and obtain a (qualitative) description 
of the global dynamics of these models. Consequently, such an analysis allows for an efficient preliminary investigation 
of these theories, suggesting what kind of models deserve further investigation. Of particular importance are those 
theories that admit solutions that have an expansion history similar to the standard ACDM model and are therefore 
worth considering as background models for a description of the growth of structure in HTG ^2^. 

In this paper, using the Dynamical Systems Approach (DSA) approach suggested by Collins and then by Ellis and 
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Wainwright (see [2l| for a wide class of cosmological models in the GR context), we develop a completely general 
scheme, which in principle allows one to analyze every fourth order gravity Lagrangian. Our study generalizes p^ . 
which considered a generic power law function of the Ricci scalar f{R) — i?" and extends the general approach given 
in a recent paper [22l |. Here a general analysis was obtained using a one -parameter description of any f{R) model, 
which unfortunately turns out to be somewhat misleading. 

The aim of this paper is to illustrate the general procedure for obtaining a phase space analysis for any analytical 
f{R) Lagrangian, which is regular enough to be well defined up to the third derivative in R. After a short preliminary 
discussion about fourth order gravity, we will discuss this general procedure, giving particular attention to clarifying 
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the differences between our approach and the one worked out in [22^ . In order to illustrate these differences and the 
problems that exist in [23], we will apply our method to two different families of Lagrangian Rp expqR and R + xR^- 
The last part of the paper is devoted to discussion and conclusions. Unless otherwise specified, we will use natural 
units (h = c — ks = SnG = 1) and the (+, — , — , — ) signature. 

II. FOURTH ORDER GRAVITY MODELS 

If one relaxes the assumption of linearity of the gravitational action the most general fourth order Lagrangian in 
an homogeneous and isotropic spacetime can be written as : 

L = V^[f{R) + CM] . (1) 
By varying equation ([1]), we obtain the fourth order field equations 

f'{R)R^, - \f{R)9f,. = f'{Rr^ (5.^5/3. - Oo^pg^.) + ff^i , (2) 

where f *5 = ^'^(v^-^m) 

prime denotes the derivative with respect to R. Standard Einstein equations 



-g og^,^ 

are immediately recovered if f{R) — R- When f'{R) ^ the equation ([2|) can be recast in the form 

G = R^i/ — -^g/luR — = + , (3) 

where 

Tjl = jj^ [f{R) - Rf'im + fiRY'^^igc^gp. - 5a/3g^.)| , (4) 

represent the stress energy tensor of an effective fluid sometimes referred to as the "curvature fluid" and 

rpM ^ ^ rpM /r\ 

represents an effective stress-energy tensor associated with standard matter. 

The conservation properties of these effective fluids are given in [lO, [1^ but it is important to stress that even if 
the effective tensor associated with the matter is not conserved, standard matter still follows the usual conservation 
equations f/^f,'" = 0. 

Let us now consider the Friedmann-Lemaitre-Robertson- Walker (FLRW) metric: 



ds"^ =dt^ -a^(t) 



' r [do + sm 'J-'-^^' 



(6) 



1 — kr^ 

For this metric the action the field equations (jH) reduce to 

^' + ^= 3F{M/'^-/]-3^/' + Mm} , (7) 
2ij + + ^ = _^ 1 1 [fR _ /] + f> _ 3ij// + , (8) 

and 
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R = -6 \^2H-' + H+-^j , (9) 

where H = d/a, f = and the "dot" is the derivative with respect to t. The system ((T]) is closed by the Bianchi 

identity for f^^: 

/i„ + 3H{^l„i + Pm) = , (10) 
which corresponds to the energy conservation equation for standard matter. 
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III. THE DYNAMICAL SYSTEM APPROACH IN FOURTH ORDER GRAVITY THEORIES 



Following early attempts (see for example [24,]), the first extensive analysis of cosmologies based on fourth order 
gravity theory using the DSA as defined in [2l[ was given in [l^. Here the phase space of the power law model 
f{R) = x^" was investigated in great detail, exact solutions were found and their stability determined. Following 
this, several authors have applied a similar approach to other types of Lagrangians [l^, and very recently this scheme 
was generalized in [2^ . 

In this paper we give a self consistent general technique that allows us to perform a dynamical system analysis of 
any analytic fourth order theory of gravity in the case of the FLRW spacetime. 

The first step in the implementation of the DSA is the definition of the variables. Following [l^, we introduce the 
general dimensionless variables : 

^ f' R ^ / o ^ Mm ^ k 

where /j,™ represents the energy density of a perfect fluid that might be present in the model. 
The cosmological equations (O are equivalent to the autonomous system : 

^ = e{2K + 2z- + {K + y + l)x) + (-3w - 1) + 2, (12) 

^ - Ve{2v + 2K + xT + A), (13) 

ze{2K -x + 2y + '^) + exyT, (14) 

dN 

^ = ne{2K~x + 2y--iw + l), (15) 

^ - 2Ke{K + y+l), (16) 

where N = | lna| is the logarithmic time and e \H\/ H. In addition, we have the constraint equation 

1 = -K - X -y + z + n, (17) 

which can be used to reduce the dimension of the system. If one chooses to eliminate the variable associated with 
the spatial curvature, we obtain 



dx 
dN 
dy 
dN 
dz 

dN 

dn 

dN 

K — z + ri — x — J/ — 1 



£ (Az - 2x'^ + {z- 2)x - 2y) + fie {x - 3w + 1), 

ye[2n + 2{z+l) + x{T -2)], (18) 
ze {2z + 2n-3x + 2)z + xe yT, 
ne(2n- 3.T + 2z-3w- 1), 



The quantity T is defined, in analogy with [24], as 



T^^. (19) 



The expression of T in terms of the dynamical variables is the key to closing the system (j39p and allows one to perform 
the analysis of the phase space. The crucial aspect to note here is that T is a function of R only, so the problem of 
obtaining T = T{x,y,z,^) is reduced to the problem of writing R = R{x,y,z,Q). This can be achieved by noting 
that the quantity 

r^-^, (20) 
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is a function of R only and can be written as 



Solving the above equation for R allows one to write R in terms of y and z and close the system (jlSp . 

In this way, once a Lagrangian has been chosen, we can in principle write the dynamical system associated with it 
using (HH]), substituting into it the appropriate form of T = T{y,z). This procedure does however require particular 
attention. For example, there are forms of the function / for which the inversion of (|2ip is highly non trivial (e.g., 
f{R) = cosh(it!)). In addition, the function T could have a non-trivial domain, admit divergences or may not be in 
the class C^, which makes the analysis of the phase space a very delicate problem. Finally, the number m of equations 
of (HH]) is always m > 3 and this implies that fourth order gravity models can admit chaotic behaviour. While this is 
not surprising, it makes the deduction of the non-local properties of the phase space a very difficult task. 

The solutions associated with the fixed points can be found by substituting the coordinates of the fixed points into 
the system 

H = aH^ , a = -l-n, + x^- z^, (22) 
. _ 3(1 + w) 
a t 

where the subscript "i" stands for the value of a generic quantity in a fixed point. This means that for a =/= the 
general solutions can be written as 

a = ao(i-<o)'/" , (24) 

Mm = a,o{t~to) . (25) 

The expression above gives the solution for the scale factor and the evolution of the energy density for every fixed 
point in which a ^ 0. When a = Q the ((22)) reduces to H = which correspond to either a static or a de Sitter 
solution. 

The solutions obtained in this way have to be considered particular solutions of the cosmological equations which 
are found by using a specific ansatz (i.e. the fixed point condition 25:]). For this reason it is important to stress that 
only direct substitution of the results derived from this approach in the cosmological equations can ensure that the 
solution is physical (i.e. it satisfies the cosmological equations ([7]) ). This check is also useful for understanding the 
nature of the solutions themselves e.g., to calculate the value of the integration constant (s). 

Also, the fact that different fixed points correspond to the same solutions is due to the fact that at the fixed points 
the different terms in the equation combine in such a way to obtain the same evolution of the scale factor. This means 
that although two solutions are the same in terms of time dependence, the physical mechanism that realizes them can 
be different 

One difference betwen our approach and the one in [2^ is that we consider a non-zero spatial curvature k. The choice 
of including a non-zero spatial curvature k has been made with the aim of obtaining a completely general analysis of a 
fourth order cosmology from the dynamical systems point of view. In addition, since most of the observational values 
for the cosmological parameters are heavily model dependent, we chose to limit as much as possible the introduction 
of priors in the analysis. However, as we write in the footnote in section 3, the limit of flat spacelike sections {K 0) 
can be obtained in a straightforward way for our examples. In fact, each fixed point is associated with a specific 
value of the variable K (i.e. a value for k) and the stability of these points is independent of the value of K. As 
matter of fact in order to consider fixed points living on the hypersurface K = 0, one has just to exclude the fixed 
points associated with K ^ 0. In addition to that, looking at the dynamical equations one realizes that K — is 
an invariant submanifold, i.e., an orbit with initial condition = will not escape the subspace K = and orbits 
with initial condition K can approach the hyperplane K — only asymptotically. As a consequence, one does 
not need to have any other information on the rest of the phase space to characterize the evolution of the orbits in 
the submanifold K — 0. The authors of [l^l proposed that the function m{r) = T(r)~^ could be used as a parameter 
associated with the choice of f{R), thus obtaining a "one parameter approach" to the dynamical systems analysis of 
/(i?) gravity. Unfortunately their method has several problems that lead to incorrect results. These problems can be 
avoided only if one considers the framewo rk p resented above. 

Let us look at this issue in more detail [28j . In the system equivalent to ([39)l is associated with the relation 

^=ril + Hr)+r)^, (26) 

which is clearly a combination of the equations for z and y. In order to ensure that the variable r and consequently 
the parameter m is constant they require the RHS of the above equation to be zero. Their solution to this problem 
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is the condition 1 + m(r) + r = 0, which is an equation for r when the function m(r) has been substituted for and is 
also the bases of their method of analysis. 

The problem here is that this equation has not been fully expressed in terms of the dynamical system variables. In 
fact, one can rewrite in the form: 

dr _ r(l + m{r) + r) 

dN ~ ^) ^ ' ^ ' 

dv 

which means that the condition — — = in fact corresponds to 

dN 

r(l +m(r) +r) 

-r\ x = (28) 

rather than 1 + m(r) + r = 0. Equation (|28p has a solution if 

a; = 0, (29) 
r = 0, (30) 

^i±^44^-0. (31) 
m[r) 

and this leads to solutions for r which are in general different from the values of r obtained from 1 + m(r) + r — Q. 
This inconsistency has major consequences for the rest of the analysis in p^ . leading to changes in the number of 
fixed points as well as their stability (see below for details). 

In fact, a more careful analysis reveals that for some of the fixed points (e.g. Pi, ...P4) the values of r obtained from 
the relation r = —y/z either cannot be determined unambiguously or do not solve the condition 1 + m(r) + r = 0, 
which is claimed to come from (|26p in p^ . 

This is a clear indication that the approach used in is both incomplete and leads to wrong conclusions. It is also 
interesting to stress that if one substitutes the expression for m in terms of the dynamical system variables in (26-29) 
of l22l . the results match the one obtained in our formalism. This implies that the reason the method described in 
[2^ fails has its roots in the attempt to describe the phase space of a whole class of fourth order theories of gravity 
with only one parameter. 

In the following we will present a number of examples of f{R) theories that can be analyzed with this method and 
we compare the results obtained with those given in [22.] . 



IV. EXAMPLES OF f{R) - LAGRANGIANS 

In this section we will show, with the help of some examples, how the DSA developed above can be applied. In 
particular we will consider the cases f{R) = exp{qR) and f{R) — R + X-R"- Since the aim of the paper is to 
provide only the general setting with which to develop the dynamical system approach in the framework of fourth 
order gravity, we will not give a detailed analysis of these models. Istead, we will limit ourselves to the finite fixed 
points, their stability and the solutions associated with them. A comparison with the results of [l^l will also be 
presented. 



A. The f{R) = R" exp{qR) case 



Let us consider the Lagrangian f{R) = RPexp{qR). As explained in the previous section, the dynamical system 
equations for this Lagrangian can be obtained by calculating the form of the parameter T. We have 



T(y,^) 



y z 



y"^ — p z"^ 



Substituting this function into ([39]) we obtain 

e [4z - + (^ - 2)a; - 2y] + 17e (x - 3w + 1), 



dx 
dN 
dy 
dN 



ye 



2n + 2z + 2- 



y"^ ~ p z"^ 



- 2x 



(32) 

(33) 
(34) 
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dz 
dN 



ze 



2z + 2n-3x + 2 • 



— p z^ 



(35) 



^ = fie {2n - 3x + 2z - 3w - 1), (36) 
K = z + n-x-y-1. (37) 

The most striking feature of this system is the fact that two of the equations have a singularity in the hypersurface 
y'^ = p z^ . This, together with the existence of the invariant submanifolds y = and z = Q heavily constrains the 
dynamics of the system. In particular, it implies that no global attractor is present, thus no general conclusion can be 
made on the behavior of the orbits without first providing information about the initial conditions. The finite fixed 
points can be obtained by setting the LHS of to zero and solving for (x, y, z, 17), the results are shown in TableU 

The solutions corresponding to these fixed points can be obtained by substituting the coordinates into the system 
([22ll and are shown in Table [2^ . The stability of the finite fixed points can be found using the Hartman-Grobman 
theorem [27]. The results are shown in Table Hill Note that some of the eigenvalues diverge for p = 0, 1. This happens 
because in the operations involved in the derivation of the stability terms p—\ and/or p appear in the denominators. 
However this is not a real pathology of the method but rather a consequence of the fact that for these two values 
of the parameter the cosmological equations assume a special form. In fact it is easy to prove that if one starts the 
calculations using these critical values of p one ends up with eigenvalues that present no divergence (23| . 

Let us now compare our results with the ones in }22|. The number of fixed points obtained for this Lagrangian, 
when if = 0, matches the ones obtained in j22]. This result can be explained by the fact that the solutions of the 
constraint equation for m (|26p coincide with the ones coming from the correct constraint equation (|27|) (the matching 
between the two systems can be obtained setting w = in Table [T|. However, when one calculates the stability of 
these points our results are strikingly different to those presented in (2^ . For example, in our general formalism it 
turns out that the fixed point M (corresponding to P5 of [22]) is a saddle for any value of the parameter p and, as 
consequence, it can represent only a transient phase in the evolution of this class of models. Instead, in [2^ the authors 
find that this point can be stable (not necessarily always a spiral) and argue that this fact prevents the existence of 
cosmic histories in which a decelerated expansion is followed by an accelerated one. From this they also conclude 
that an entire subclass of these models (m — m{p) > 0) can be ruled out. Our results show clearly that this is not 
the case. Another example is the point A4 corresponding to Pg of [13] . In the authors find that this point can 
be stable or a saddle as we do, but the intervals of values of the parameters for which this happens are different (see 
Table Uni. As explained above, the reason behind these differences is the fact that the method used in leads to 
incorrect results when, like in this case, there is no unambiguous way of determining the parameter r = —y/z from 
the coordinates of the fixed points. Consequently the conclusions in [24] relating to the properties of these points are 
incorrect and have no physical meaning. 



TABLE L Fixed points of i?'' exp(gi?).The superscript represents a point corresponding to a double solution. 
Point Coordinates {x, y, z, Q) K 



A 


(0,0,0,0) 


-1 


B 


(-1,0,0,0) 





C 


(-l-3u;,0, 0, -l-Sw) 


-1 


V 


(l-3w,0,0, 2-3w) 
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(2,0,2,0) 


-1 




(1,-2,0,0) 





g 


(0,-2,-1,0) 





H 


(4,0,5,0) 





I* 


(-3(1 + u>), -2, 0,-4-310) 





L 


(2-2p,2p(l-p),2-2p,0) 


2p(p - 1) - 1 


M 


/ 4-2p (5-4p)p 5-4p p,\ 
I l-2p' 2p^-3p+l' (p-l)(2p-l) ' '-'y 





N 


/ -3{1+jd){p-1) 3(l+m)-4p -4p+3iu+3 p(9iu-2p{3ui+4) + 13)-3('uj+l) \ 





I p ' 2p ' 2p^ ' 2p^ J 
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TABLE II: Solutions associated with the fixed points of R''exp{qR). The solutions are physical only in the intervals of p 
mentioned in the last column. 



Point 



A 
B 
C 
V 
£ 

J"* 
Q 

n 



M 



Scale Factor 



Energy Density 



Physical 



a{t) = {t- to) 
a{t) = ao{t — to)^^"^ 

a{t) = (t - to) 
a{t) = ao{t — toY^^ 

a{t) = {t- to) 
a{t) — ao, 



{t - to) 



6v^ (^- 

ao {t-tof^ 



to) 



a{t) — ao exp 
a(t) = ao, 
a{t) — ao cxp 

ait) 
a{t) = ao, 
a{t) = ao exp 
a{t) = {t- to) 1/1 - 2p(p - 1) 

a{t) = ao(t- to) 2-p 
a{t) = ao (t-to)^^^ 








3(2p2 



-3p+ll(,<, + l) 



% = Hmojt — to)~^'' 



P>1 

p>2 
p>l 
p>2 
p>l 

p>0 
3,y ^ Vp > 3, 

p > 

p< |,g>OVp> f,g<0 
p> 2 
p> 



p< |,g>OVp> 2 g<0 



p< |,g>OVp> 



1<P<I+ 

1 



3 ' 

V3 
2 ' 2 



2 g < 



3(tii+l) 



-12 
2 ' ' 4 



0) 



TABLE III: The stability associated with the fixed points in the model B'' cxpiqR). With the index we have indicated the 
attractive nature of the spiral points. 

Point Stability 



A saddle 



g r repellor < w < 2/3 

[ saddle otherwise 
C saddle 

^ ( repellor 2/3 < w < 1 
1^ saddle otherwise 



£ saddle 
saddle 

attractor 0<w;<lU2<p<|| 
spiral+ 0<w<luf|<p<4 
saddle otherwise 
7i saddle 
I non hyperbolic 

( attractor i-^^<p<OV|<p<i + ^ 
<^ spiral+ < p < | 

saddle otherwise 



M 



saaaie otnerwise 
attractor p < i(l - a/3) V §(1 + a/3) < p < 2 
saddle otherwise 



J\f saddle 
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B. The case f{R) ^ R + xBJ" 

Let us discuss now the case of a Lagrangian corresponding to a power law correction of the Hilbert - Einstein gravity 
Lagrangian f{R) = i? + X-R"- In this case, the characteristic function T(y, z) reads : 

T(y,^) = ^ (38) 
n{z - y) 

and substituting this relation into the system of equations (pS)) one obtains 



dx 

In 



2x-' + {z-2)x-2y + 4:Z + n{x-3w + l), (39) 



^ = ye [2n + 2{z + 1) + - 2x], (40) 

aN n(z — y) 

% = 7z£(2z + 2r>-3x + 2) + £ (41) 
aiV n(z — y) 



dfl 

In 

K = z + ^-x-y-l. (43) 



f7£(2r2-3j: + 2z-3w-l), (42) 



As in the case of /(i?) = i?^exp(gi?), the system is divergent on a hypersurface (this time y = z) but it admits only 
one invariant submanifold, namely y = 0. This, again, implies that no global attractor is present and no general 
conclusion can be made on the behavior of the orbits without giving information about the initial conditions. The 
finite fixed points, their stability and the solutions corresponding to them are summarized in Tables HVl IVl and IVII 

As before our results are different from those given in [24I . First of all, our set of fixed points do not coincide with 
the ones presented in [2^ . In particular, in our analysis there is no fixed point corresponding to Pc,a. Again, the 
reason for this difference is to be found in the constraint equation (j26p . which in this case gives the incorrect set of 
solutions and therefore affects the set of fixed points. In fact, if one substitutes the expression for m{r) of [1^ in 
terms of the coordinates in equations (34)- (39), it is easy to verify that two of these equations diverge at this point. 

The differences between the results in our approach and the one presented in [23| are even more evident when the 
stability analysis is considered. For example, the point £, corresponding to Pi, is always a saddle, except into the 
region < n < 2 when it is attractive. This behavior is recovered in only for —2 < n < —41/25. Also, points G 
(corresponding to P4 of [l^l) and V (corresponding to P3 of [l^]), which in our approach are always saddles in the 
dust case, are always repellers in Finally, also the stability of I corresponding to Pg appears to be different from 
the one presented in [22l |. 



TABLE IV; Coordinate of the finite fixed points for R + xR" gravity. 

Point Coordinates {x,y,z,Q,) K 

A (0,0,0,0) -1 

B (-1,0,0,0) 

C (-1 - 3«;, 0,0,-1 - 3u>) -1 

V (l-Sw, 0,0,2-310) 

S (0,-2,-1,0) 

T (2,0,2,0) -1 

G (4,0,5,0) 

U (2(1 - n),2n(n - 1),2(1 - 7i),0) 2n(n - 1) - 1 

q- f 2(n-2) (5-4n)n 5-4n r.\ rj 

V 2n-l ' 2n^-3n+l' 2n^-3n+l'"/ ^ 
f ( 3(n-l)(ro + l) -4n + 3u) + 3 -4n+3ro+3 -2(3u)+4)n'' +(9m + 13)n-3(m + l) \ 

' n ' 2n ' 2n^ ' 2n^ ' ^ 
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TABLE V: The stability of the fixed points in the model R + x^"- The quantities Bi related to the fixed point C, represent 
some non fractional numerical values {B\ ~ 1.220, B\ ~ 1.224, B3 ~ 1.470). 



Point Stability 



A 
B 
C 
V 



T 
Q 

n 



C 



saddle 

repellor < ui < 2/3 

saddle otherwise 
saddle 

J repellor 2/3 < < 1 

1 saddle otherwise 



attractor 
spiral^ 
saddle 
saddle 
saddle 

attractor 
spiral^ 
saddle 
attractor 

repeller 

saddle 
tj; = 0,1/3 

w = 1 



0<n< § 
otherwise 



^(1 



- ^3) < n < 
< n < 1 
otherwise 
n < i(l - y3)Un > 2, 
l<n< = 0,1/3), 

l<n< i(ll + ^),(w = l) 
otherwise, 
saddle, 

repellor Bi < n < 52 U B3 < n < | , 



saddle 



otherwise 



TABLE VI: Solutions associated to the fixed points of R + xR"- The solutions are physical only in the intervals of p mentioned 
in the last column. 



Point 



A 
B 
C 
V 

£* 

T 
Q 
H 

X* 

c 



Scale Factor 



a{t) = {t- to) 
a{t) = ao(t- to)^^^ 

a{t) = (t - to) 
a{t) = ao(t- to)^^^ 

a{t) = ao, 

a(t) =aoexp[±2^/3x^(2-3n)^(^-^o)] 

~ 2(l-n) 

a{t) = {t- to) 
a{t) = ao{t~ to)^^^ 
a(t) = ^1 - 2n{n - 1) {t - to) 

2ti^-3ti + 1 

a{t) ^ ao {t - to) 2-™ 
a{t) ^ao {t-to)^H^ 



Energy Density 



Physical 



3(2n2-3»i+l)(ii. + l) 
/^m — f^m oi " 

Mm ~ tJ-m,o{t — tp)^^ 



n > 1 
n > 1 
n > 1 
n > 1 
n > 
n< |,x > V 
n > |,x < 
n > 1 
n > 1 



l<n> i+ , 



n = 2 ) Mm ,0 = 

non physical 



V3 





CONCLUSIONS 



In this paper we have presented a general formalism that allows one to apply DSA to a generic fourth order 
Lagrangian. The crucial point of this method is to express the two characteristic functions [2^ : 



T = 



f 

Rf" ' 



Rr_ 



(44) 
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in terms of the dynamical variables, which, in principle, allows one to obtain a closed autonomous system for any 
Lagrangian density /(i?). 

The resulting general system admits many interesting features, but is very difficult to analyze without specifying 
the function T (i.e. the form of f{R)). Consequently, a "one parameter" approach can lead to a number of misleading 
results. 

Even after substituting for T, the dynamical system analysis is still very delicate; in fact, T could be discontinuous, 
admit singularities or generate additional invariant submanifolds that influence deeply the stability of the fixed points 
as well as the global evolution of the orbits. 

After describing the method, we applied it to two classes of fourth order gravity models: R + xi?" and i?^ exp{qR), 
finding some very interesting preliminary results for the finite phase space. Both these models have fixed points with 
corresponding solutions that admit accelerated expansion and, consequently can model either infiation or dark energy 
eras (or both). In addition, there are other fixed points which are linked to phases of decelerated expansion which can 
in principle allow for structure formation. These latter solutions are not physical for every value of their parameters, 
but this is not necessarily a problem. In fact, in order to obtain a Friedmann cosmology evolving towards a dark 
energy era, these points are required to be unstable i.e., cosmic histories coast past them for a period which depends 
on the initial conditions. This means that the general integral of the cosmological equations corresponding to such 
an orbit will only approximate the fixed point solution and this approximate behavior might still allow structures to 
form. 

It is also important to mention the fact that even if one has the desired fixed points and desired stability, this does 
not necessarily imply that there is an orbit connecting them. This is due to the presence of singular and invariant 
submanifolds that effectively divide the phase space into independent sectors. Of course one can implement further 
constraints on the parameters in order to have all the interesting points in a single connected sector, but this is still 
not sufficient to guarantee that an orbit would connect them. The situation is made worse by the fact that, since the 
phase space is of dimension higher than three, chaotic behavior can also occur. It is clear then, that any statement 
on the global behavior of the orbits is only reliable if an accurate numerical analysis is performed. However, these 
issues (and others) will be investigated in more detail in a series of forthcoming papers. 

A final comment is needed regarding the differences between our results and the ones given in (2^ . Even if the 
introduction of T and r, was suggested for the first time in that paper, the results above (and in particular the 
existence of a viable matter era) are in disagreement with the ones given in that paper. The reason is that the 
authors of ^] used "a one parameter description" in order to deal with ([55]) in general. We were able to prove that, 
unfortunately, not only are the equations given in 22] incomplete, but also that the method also gives both incorrect 
and misleading conclusions. 



References 



[1] V. Sahni and A. A. Starobinsky, Int. J. Mod. Phys. D 9, 373 (2000); T. Padmanabhan, Phys. Kept. 380, 235 (2003); P. 

J. E. Peebles, B. Ratra, Rev. Mod. Phys. 75, 559 (2003); V. Sahni, Lect. Notes Phys. 653, 141 (2004); E. J. Copeland, 

M. Sami, S. Tsujikawa, arXiv:hep-th/0603057 
[2] S. Perlmutter et al, Astrophys. J. 517, 565 (1999); A. G. Riess et al, Astron. J. 116, 1009 (1998); J. L. Tonry et al, 

Astrophys. J. 594, 1 (2003); R. A. Knop et ai, Astrophys. J. 598, 102 (2003); A. G. Riess et al. Astrophys. J. 607, 665 

(2004); S. Perlmutter et al. Astrophys. J. 517, 565 (1999); Astron. Astrophys. 447, 31 (2006). 
[3] D. N. Spergel et al. Astrophys. J. Suppl. 148, 175 (2003); D. N. Spergel et al. arXiv: astro- ph/0603449 
[4] M. Tegmark et al., Phys. Rev. D 69, 103501 (2004); U. Seljak et al., Phys. Rev. D 71, 103515 (2005); S. Cole et al, Mon. 

Not. Roy Astron. Soc. 362, 505 (2005). 
[5] D. J. Eisenstein et al., Astrophys. J. 633, 560 (2005); C. Blake, D. Parkinson, B. Bassett, K. Glazebrook, M. Kunz and 

R. C. Nichol, Mon. Not. Roy Astron. Soc. 365, 255 (2006). 
[6] B. Jain, A. Taylor, Phys. Rev. Lett. 91, 141302 (2003). 
[7] P. Astier et al, Astron. Astrophys. 447, 31 (2006). 

[8] S. Capozziello, V. F. Cardone, S. Carloni, A. Troisi, Int. J. Mod. Phys. D 12, 1969 (2003). 
[9] A. A. Starobinsky, Phys. Lett. B 91, 99 (1980). 
[10] S. Capozziello, S. Carloni and A. Troisi, Recent Res. Dev el. Astronomy & Astrophysics 1, 625 (2003), arXiv; 
astro-ph/0303041 

[11] S. Capozziello, V. F. Cardone, A. Troisi, JCAP 0608, 001 (2006). 



11 



[12] K. i. Maeda and N. Ohta, Phys. Lett. B 597, 400 (2004) arXiv:hep-th/0405205 , K. i. Maeda and N. Ohta, Phys. Rev. D 
71, 063520 (2005) arXiv:hep-th/0411093 , N. Ohta, Int. J. Mod. Phys. A 20, 1 (2005) arXiv:hep-th/0411230; , K. Akune, 



[14 
[15 

[16: 
[17: 
[is: 

[19 
[20 

[21 



[22' 
[23^ 
[24 
[25^ 
[26 
[27 
[28 



K. i. Maeda and N. Ohta, Phy s. Rev. D 73, 103506 (2006) [arXiv:hep-th /0602242 

[13] S. Nojiri and S. D. Odintsov, arXiv:0810.1557 [hep-th]; S. Nojiri and S. D. Odintsov, 'arXiv:0807.0685' [hep-th]; S. No- 
jiri and S. D. Odintsov, arXiv:0801.4843 [astro-ph]; G. Cognola, E. Ehzalde, S. Nojiri, S. D. Odintsov, L. Sebastian! 
and S. Zerbini, Phys. Rev. D 77, 046009 (2008) arXiv:07rO017 [hep-th]]; S. Nojiri and S. D. Odintsov, J. Phys. 
Conf. Ser. 66, 012005 (2007) arXiv:hep-th/0611071 ; S. Nojiri and S. D. Odintsov, Phys. Rev. D 74, 086005 (2006) 
[arXiv!hep- th /0608008 ; S. Nojiri and S. D. Odintsov, eConf C0602061 (2006) 06 [Int. J. Geom. Meth. Mod. Phys. 4 
(2007) 115] arXiv:hep-th/0601213 ; S. Nojiri and S. D. Odintsov, Phys Lett. B 631 , 1 (2005) arXiv:hep-th/0508049 ; 
S. Capozziello, S. Nojiri and S. D. Odintsov, Phys. Lett. B 632, 597 (2006) [arXiv:hep-t"h /0507182 ; S. Nojiri, S. D. Odintsov 
and M. Sasaki, Phys. Rev. D 71, 123509 (2005) arXiv:hep-th/0504052 ;M. C. B. Abdalla, S. Nojiri and S. D. Odintsov, 
Class. Quant. Grav. 22, L35 (2005) arXiv:hep-th/0409177 ; S. Nojiri and S. D. Odintsov, Gen. Rel. Grav. 36, 1765 (2004) 
arXiv:hep-th/0308176 ; S. Nojiri, S. D. Odintsov, Phys. Lett. B 576, 5, (2003); S. Nojiri, S. D. Odintsov, Phys. Rev. D 
68, 12352 (2003); 

S. Capozziello, V. F. Cardone, A. Troisi, Mon. Not. Roy. Astron. Soc. 375, 1423 (2007). 

R. Kerner, Gen. Relativ. Grav. 14 453 (1982) ; J. P. Duruisseau, R. Kerner, Glass. Quantum Grav. 3 817 (1986). 
P. Teyssandier, Glass. Quantum Grav. 6 219 (1989). 

G. Magnano, M. Ferraris and M. Francavigha Gen. Relativ. Grav. 19 465 (1987). 
S. Carloni, P. Dunsby, S. Capozziello S, Troisi A 2005 Class. Quantum Grav. 22, 4839. 
S Carloni, J. Leach S. Capozziello, P. Dunsby Class. Quantum Grav. 25 035008 (2008). 

K. N . Ananda, S. Carloni and P. K. S. Dunsby, arXiv:0809.367 3 [astro-ph]; K. N. Ananda, S. Carloni and P. K. S. Dunsby, 
[axXiV :0812.2028 [astro-ph] 

Dynamical System in Cosmology edited by J. Wainwright and G F R Ellis (Cambridge: Cambridge Univ. Press 1997) and 
references therein 

L. Amendola, R. Gannouji, D. Polarski and S. Tsujikawa, Phys. Rev. D 75 083504 (2007). 
S. Carloni, A. Troisi and P. K. S. Dunsby, arXiv:0706.0452 [gr-qc] (2007). 
S. Capozziello, F. Occhionero and L. Amendola, Int. J. Mod. Phys. D 1, 615 (1993). 
S. Carloni, K. N. Ananda, P. K. S. Dunsby and M. E. S. Abdelwahab. .,arXiv:0812.221 1 [astro-ph]. 
J. D. Barrow and S. Hervik, Phys. Rev. D 74, 124017 (2006). 
P. Hartmann Ordinary differential equations (New York Wiley, 1964) 

It is important to note that in 22] the signature is not the same of the one used here (e.g -,+,+,+ instead of 4-,-,-,-) and 



the definition of the variables are slightly different. The transformation from one variable to another is as follows; 

X — + — xi, y — > — X3, z — * X2, K ^ 0, w ^ 0. 



However, as expected, this does not affect our conclusions. 
[29] Note that even if the parameter q is not present in the dynamical equations it appears in the solutions because we have 
calculated the integration constants via direct substitution in the cosmological equations. 



